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We show that a large number of ions stored in a Penning trap, and forming a 2D Coulomb crystal, 
provides an almost ideal system for scalable quantum computation and quantum simulation. In 
particular, the coupling of the internal states to the motion of the ions transverse to the crystal 
plane, allows one to implement two qubit quantum gates. We analyze in detail the decoherence 
induced by anharmonic couplings with in-plane hot vibrational modes, and show that very high 
gate fidelities can be achieved with current experimental set-ups. 



The search for a physical system where quantum com- 
putation is feasible is at the focus of an intense theo- 
retical and experimental activity [1]. Ion traps are by 
now among the most promising candidates for a many- 
qubit quantum processor. In this system, qubits are 
stored in internal electronic states, and collective vi- 
brational modes of the ions allow us to induce quan- 
tum gates between them [2]. Following this idea, the 
building blocks for quantum computation have already 
been demonstrated in experiments with a few qubits [3]. 
Most of the current efforts to scale up the size of cur- 
rent ion quantum processors, rely on the fabrication of 
arrays of microtraps [4] , in which a large number of ions 
can be stored and shuttled. Even though an astonishing 
progress has been achieved in this direction in the last 
years, the scalability of this system still demands techni- 
cal advances in microfabrication and trap design [5]. 

Penning traps provide us with an alternative trapping 
scheme, where a large number of ions (10 4 - 10 6 ) can 
be confined by a potential with approximate cylindrical 
symmetry [6]. Axial confinement is induced by a static 
electric field, whereas radial confinement is a result of 
the rotation of the ions under an axial magnetic field. 
If the axial confinement is strong enough, ions arrange 
themselves in a triangular lattice on a single plane, which 
corresponds to a classical two dimensional (2D) Wigner 
crystal. The appeal of this system lies on the fact that 
ions are naturally ordered in a 2D regular array, with- 
out the need of individual micropotentials. Furthermore, 
ions are separated by distances of the order of tens of mi- 
crons, such that they are individually addressable by op- 
tical means [7]. Thus, ions in Penning traps may appear 
as ideally suited for quantum computation and quantum 
simulation. However, this system has never been con- 
sidered for this task [8]. First, because the complicated 
vibrational level structure of the 2D crystal makes it dif- 
ficult to apply here schemes that require resolution of 
single vibrational modes. Beside that, typical schemes 
usually rely on the coupling of qubits to modes in di- 
rections parallel to the crystal. In current experiments 
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FIG. 1: (color online). Quantum gate in a 2D Coulomb crys- 
tal: standing-waves induce a state-dependent dipole force on 
two nearest neighbors in a triangular lattice. 



with Penning traps, Doppler cooling of ions has reached 
temperatures of at most 1 mK, which implies occupation 
numbers of 10 2 - 10 3 in the in-plane vibrational modes, 
so that, it seems not to be possible to use them for quan- 
tum operations. 

In this Letter we show how to circumvent these prob- 
lems by exploiting the ions' motion along the axial direc- 
tion (Fig. 1). This approach benefits from the high ax- 
ial confinement frequencies (and thus smaller occupation 
numbers at finite temperature), as well as from the fact 
that ions are weakly coupled in this direction, something 
that enormously simplifies the description of the ions' 
motion. In particular, we show how it is possible to carry 
out two-qubit gates between ions with high fidelities by 
performing a careful analysis of the main sources of de- 
coherence. We emphasize that the results derived here 
also imply that this system is ideally suited for quantum 
simulations, which may be specially interesting due to 
the fact that ions are displayed in a triangular structure 
and favor the simulation of magnetic frustrated systems. 

The main source of decoherence in our scheme is due to 
the anharmonic terms in the Coulomb interaction, which 
induce a coupling between axial motion and in-plane hot 
vibrational modes, and lead to a residual qubit-phonon 
coupling. The description of such decoherence poses an 
involved theoretical problem, because of the large num- 
ber of vibrational modes that participate in the process. 
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However, it gets simplified due to the fact that the envi- 
ronment is in a gaussian state, which allows us to sim- 
plify the calculation of correlation functions that appear 
in finite-temperature, time-dependent perturbation the- 
ory. We show that in a range of parameters where ax- 
ial confinement is large enough, the error induced in the 
quantum gate is very small. Furthermore, the adjust- 
ment of the gate time allows us to correct for the influ- 
ence of the phonon environment, and decrease the gate 
error by more than one order of magnitude. 

Let us consider a system of N ions forming a 2D 
Coulomb crystal in a Penning trap. We study, for con- 
creteness, the performance of a "pushing gate" [9] in the 
axial direction, an approach that has the advantage that 
single vibrational modes do not have to be resolved, and 
the gate can operate at finite temperature. Since we are 
interested in estimating the consequences of decoherence, 
we can neglect finite size effects, and describe the crys- 
tal by a regular triangular lattice with periodic boundary 
conditions. Let us assign the z direction to the axis of 
the trap, such that ions occupy equilibrium positions in 
the x-y plane: 



R" = (riai+r 2 a 2 )c?o, r,- = 1, 



(1) 



where ai = (1,0,0), a 2 = (1/2,^3/2,0), and d is the 
distance between ions. The potential is given by a trap- 
ping term, plus the Coulomb repulsion: 



V = Vrrap + Vcoul) 



Vi 



Coul 



= \Y. 



R2 Rg + R r R s 



(2) 



R r = (X r , Y T ,Z r ), are the ions' coordinates with respect 
to the equilibrium positions, and Vrr ap is the harmonic 
trapping potential, with frequencies uij in each spatial 
direction, j = x,y,z. Note that in a Penning trap, (2) 
corresponds to the potential in a frame rotating with the 
ion crystal [10]. 

In the harmonic approximation, Vcoul is expanded up 
to second order in R r , and the axial (z) and in-plane 
(x,y) modes are independent. The vibrational harmonic 
Hamiltonian is given by: 



H. 



(0) 



vib 



(3) 



q is the phonon wave-vector, and A runs over the three 
possible polarizations: A = z (axial modes), and the two 
in-plane modes (A = ||,_L) corresponding to longitudi- 
nal and transverse modes. Local displacements can be 
expressed in terms of collective coordinates: 



R r = (l/ViV)]>>Xe lqr , 

A,q 



(4) 



where eq are the polarization vectors, which are eigen- 
statcs of the Fourier transform of the harmonic ion-ion 



interaction [11]: 
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For each q, are given by the eigenvectors of the 3x3 
matrix -^/f2 q . 

The "pushing gate" works by coupling the internal 
states of two ions, which we label 1 and 2, to the ax- 
ial (z) motion [9, 12]. An off-resonant standing-wave 
induces a force which displaces the position of ions 1 and 
2, in a direction which depends on their internal state 
[13] (Fig. 1): 



(6) 



i=i,2 



where <rj are operators acting on the internal states of 
ions at sites rj. To simplify our calculations, we consider 
the following time-dependent force: F(t) = Fe^ lT ^^ 2 
(— oo < t < oo). However, our results arc qualitatively 
valid for other pulse shapes with gate time 1 /T. Trapping 
parameters are chosen such that Z = e 2 jmuj 2 ^ <C 1 
[13]. In this limit, ions moving in the axial direction 
can be considered as independent harmonic oscillators 
weakly coupled by the Coulomb interaction, akin to the 
case of microtrap arrays. Furthermore, we will be in the 
adiabatic limit, defined here by: 

£ z = 8(T/cu z ) 2 (FZ Q /hu Jz ) 2 (2n z + 1) « 1, (7) 

where Z is the axial ground state size, and n z is the 
mean axial phonon number. £ z is, indeed, the error in- 
duced in the quantum gate by nonadiabatic effects in 
the axial vibrational degrees of freedom [12]. Under con- 
dition (7), internal states end up being decoupled from 
the axial motion after the gate, and follow the unitary 
evolution given by: 



1/ = e -*.M*J(*K»S } J(t) = 20 x 



F(t)Zo 

hw z 



(8) 



In the following we consider the choice T = J(0)w/8, 
such that U g corresponds, up to local operations, to a 
sign gate. 

Unfortunately, the axial motion is coupled to the in- 
plane modes by anharmonic terms in Vcoul, being the 
lowest order ones of the form XZ 2 , X 2 Z 2 , that is, 
quadratic in the axial coordinates. Since resonances be- 
tween axial and in-plane vibrational frequencies lead to 
divergences in the correlation functions of these terms, 
the effect of anharmonicities is reduced if 3> 2 Wq 1 . 
The axial vibrational bandwidth is proportional to f3 z , 
and, thus, this condition can be imposed by choosing a 
thight enough axial confinement (sec [14]). In the case 
of N — 10 4 ions, uj z w 50 ui xy is enough to ensure this 
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FIG. 2: Spectrum of the vibrational modes with wavevector 
q along the x direction. N = 10 4 ions, ui z = 50w xy , and we 
consider periodic boundary conditions. 

condition (see Fig. 2), so that we will assume this ratio 
in all the examples presented along this work. 

Let us now quantify the loss of fidelity induced by an- 
harmonic couplings. The force (6) displaces the ions in 
the axial direction by Z(t) — 2F(t)ZQ /fiw z . In the limit 
(7), ions 1 and 2 follow adiabatically the displacement 
induced by the force, so that one can neglect the fluctu- 
ations in the coordinates Z rj , and replace them by their 
ground state average, Z(t)a Z j, when computing anhar- 
monic corrections. The anharmonic energy dependence 
on the ions' position results in the following coupling be- 
tween internal states and in-plane modes: 

H$(t) = Z{tf (A T +R^ 2 B R li2 ) . (9) 

Ri ; 2 (R? 2) i s the vector given by the in-plane compo- 
nents of R ri — R r2 (R" — R" 2 ) . A and B are third and 
fourth order anharmonic Coulomb interaction terms, re- 
spectively: 

A-3^R M , B-^l-- ^ J. (10) 

Note that i? xy c excites the in-plane phonons depend- 
ing on the ions' internal states, and, thus, entangles the 
qubits with the environment. 

Assume that the internal states are initially in a pure 
state: |tf> = £ Q c Q |a) (|a> = |00), |01), |10), |11)), such 
that the initial density matrix of the total system is given 
by: Pi = I^X^I ® Pxy, W1 th Pxy a thermal phonon state. 
After the action of the quantum gate, the final density 
matrix is: 

PS = E c «4 {Uz\a)(J)\Ul) (t/ x >° y «)t) , (11) 

where £/ xy , is the evolution operator of the in-plane 
phonons, f7 xy , projected in the internal state a. In (11) 



we have assumed that the axial modes can be described 
classically in -ff xy c due to the adiabaticity of the quan- 
tum gate, such that we can factorize the evolution oper- 
ator. We define the (worst-case) reduced fidelity, T, of 
the quantum gate by the overlap of the qubits' reduced 
density matrix obtained from (11), with the qubit quan- 
tum state after a perfect gate, minimized over all possi- 
ble two-qubit initial states. Note that £/ xy = £/ xy = l4y , 

and C/° y = = Ulf ] , where U^J and U^f ] , are the in 
plane phonon evolution operator in the absence and pres- 
ence of anharmonic couplings, respectively. Thus, the fi- 
delity is completely determined by the following complex 
quantity T: 

T = tr{|01)(00|p i }=tr xy {c4 a y h ]p xy C/ x O] t } 

= tr xy |Texp (-ipHgWdS) p° y |, (12) 

which corresponds to the mean value of the evolution 
operator llffl , in the interaction picture with respect to 
H^y- Note that (12) is similar to the expression found 
in the theory developed in [16]. The worst-case error is 
E = \-T = (l-5R(.F))/2. However, note that the fidelity 
can be improved, since, according to Eq. (9), the spin- 
phonon coupling depends on the operator (1 — cferf). 
For this reason, a correction of the gate time allows us 
to cancel the phase of T, and define £' = (1 — |.F|)/2, 
which quantifies the worst-case error after the correct 
calibration of the gate duration. 

Eq. (12) is a good starting point for perturbation the- 
ory, which can be carried out by expanding the time- 
ordered exponential. First of all, let us study the scal- 
ing of the terms appearing in if xy , and check whether 
a peturbative approach is indeed justified. The anhar- 
monic coupling (9) can be rewritten in terms of collective 
variables in the interaction picture: 

<(r)^F q A (r)^(r) + £ G$ (r) R^(r) R* (r) , 

A,q A.A',q,k 

(13) 

According to Eqs. (8, 9, 10), each term scales like: 

5 s^~'(£) , "' r (£)'- «"» 

where X is the size of the ground state in the radial 
trapping potential, and by O, we mean the square root 
of the variance of the operator in the ground state. Since 
the evolution of the in plane modes is governed by the 
in-plane trapping frequency, u> xy , we expect that terms in 
perturbation theory scale according to r/aj xy in Eq. (14). 
We study two sets of experimental parameters with w xy = 
20 (200) kHz, and uj z = 1 (10) MHz, which implies d a = 
46.8 (10.1) /urn. Let us consider N = 10 4 ions, such that 
both values lead to (3 Z — 3.8 • 10~ 3 . Typical temperatures 
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of T = 1 mK can be reached after Doppler cooling, such 
that mean phonon numbers are n z rs 20 (2) in the axial 
modes, and n xy w 10 3 (10 2 ), in the center of mass in- 
plane modes. If we choose (FZo/hui z ) — 0.234, such that 
£ z < 10 -5 , we get T/aj xy = 5 10~ 2 . Taking into account 
that Xo/do is small (as it should be, if the assumption of 
harmonic vibrational modes is valid) , all terms are small 
in (14). In particular, X /d = 3.6 ■ 10~ 3 , and 5.0 • 10~ 3 , 
in the cases uj xy = 20, and 200 kHz, respectively. 

To quantify precisely the error in the quantum gate, 
we proceed as follows (see the Appendix for the de- 
tails). We expand the time ordered exponential (12) up 
to fourth order in H*~, and keep all the terms up to 

order (Xo/do) 4 . Each contribution can be expressed in 
terms of time integrations of correlation functions of col- 
lective coordinates R*. Resonances appear, for example, 

in those terms of second order in G X 't which are of the 

q,K 

form (i2£(ti)i^'(ii)i£ q (*2)i£k(*2)) The problem gets 
numerically tractable by the application of Wick's the- 
orem, which allows us to express real time high order 
correlation functions in terms of the two-operator corre- 
lation function. Note that Wick's theorem can be applied 
here, even when we arc at finite temperatures, because 
of the statistical properties of the phonon thermal (gaus- 
sian) state. Finally, a summation of sets of diagrams up 
to infinite order is possible by means of the linked clus- 
ter theorem, such that (12) takes finally the form of the 
exponential of low order irreducible diagrams [17]. 

The results of our calculation are presented in Fig. 3, 
where we show our results for the two sets of experimental 
parameters discussed above. Note that at high tempera- 
tures, like those ones that occur in current experiments 
with Penning traps, we get a quadratic dependence with 
temperature, which can be explained by the dominant 
contribution of terms that are second order in the cou- 
pling G q '£ . The correction of the phase by adjustment of 
the gate time allows us to reduce the error by more than 
one order of magnitude. Note that even with the highest 
temperatures considered here, which correspond to the 
limit of Doppler cooling, anharmonic terms induce very 
small errors. With the range of parameters considered in 
this work, we get rates T = 1, 10 kHz. A limitation in 
the number of quantum gates due to heating, reduces the 
number of gates to a maximum of w 10 2 , with present 
heating rates, however, this quantity is amenable to be 
improved by increasing the quality of the vacuum in the 
trap [18] . 

Finally, our proposal could also be used for the quan- 
tum simulation of interacting spin-systems, by applying 
a state dependent force to all the ions at the same time. 
In this way, as shown in Ref. [19], an antifcrromagnetic 
Ising interaction is induced between the internal states, 
which behave like effective spins. If we add a transverse 
field of the form (0/2) J2 r cr*, by means of a global car- 
rier transition, then this experimental set-up allows us 
to study the rich phenomenology of quantum frustration 
in triangular lattices [20]. Note that our analysis of the 




FIG. 3: Continuous lines: £ (error without phase correction). 
Dashed lines: £' (error with correction of the gate time. Thick 
lines and thin lines, correspond to aj xy = 20 kHz and cj xy = 
200 kHz, respectively. N = 10 4 ions and r/oj xy = 0.05. 

decoherence induced by low-energy vibrational modes, 
also implies the viability of this approach, since the ef- 
fective spin-spin interactions are an always-on version 
of the qubit-qubit coupling induced during the quantum 
gate. 
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QUEST, SCALA, and Marie Curie under contract 
MEIF-2004-010350. 



APPENDIX A: PROPERTIES OF A WIGNER 2D 
CRYSTAL 

(In this appendix, we describe in more detail the prop- 
erties of a 2D Wigner crystal, as well as the calculation 
of the decoherence in the quantum gate due to coupling 
with low-energy vibrational modes. For the sake of read- 
ability, we repeat a few of the equations already presented 
in the main text.) 

We express the potential in the rotating frame in terms 
of the equilibrium positions and the coordinates of 
the ions with respect to the equilibrium positons, R r : 

V = Vboul + VlY a p , 

v = 1 V 

"Coul 2 ^ |R0_ R + R r _R g | ' 

V Tl . ap = im^ 2 £(Z r ° + Z r f + 

r 

+ imc r 2 ^((X r °+X r ) 2 +(y r +F r ) 2 )-(Al) 

r 

In the case uj z ^> w xy , the Coulomb crystal is a single 
plane, which corresponds to a 2D Wigner crystal. Under 
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these conditions, ions arrange themselves in the triangu- 
lar lattice generated by (1). 

The theoretical description of this system becomes 
much simpler if we consider periodic boundary condi- 
tions, something that is a good approximation, since we 
deal with large number of particles. The equilibrium po- 
sitions in a 2D Wigner crystal are given by a triangular 
lattice generated by: 



lattice). Local coordinates can be finally written like: 



e - iqr e A 



A,q 



' m 



Huh 



( a A,q + a A,-q) , 
(4,q - «A,-q) , (A7) 



The reciprocal lattice is generated by the following vec- 
tors: 



bi = (l,-l/V3,0), b 2 = (0,2/^3,0), 



(A3) 



ai = (1,0,0), &2 = (1/2, Vs/2, 0). (A2) where e A determine the direction of vibration of the 

mode, and the index A can take three values which corre- 
spond to longitudinal or transverse modes in the in-plane 
direction, and vibration perpendicular to the Coulomb 
crystal plane. 

The anharmonic corrections to the vibrational Hamil- 
tonian that are relevant for us, are those which couple the 
axial (z) to the in-plane motion (x, y), because they will 
lead eventually to the qubit-phonon couplings presented 
in Eq. (9). Let us writte them explicitly: 

= \ E4s (R r , s r(R r , s r + 

r ,s 

i— x,y 

+ \ B H (Rr,s) , (Rr,s) i (Rr,s) Z . (A8) 



The basis of the reciprocal lattice is determined by the 
relations: (a)j • (b)j — Sij. They are the allowed direc- 
tions of propagation given by the periodicity of the ion 
lattice. 

Since we are considering a lattice with regular spacing, 
we can easily determine the ions' equilibrium positions. 
In this simple case, this task is reduced to finding the 
equilibrium distance between ions, do, which we prefer 
to express in terms of the adimensional parameter /3 xy : 



/3xy(£) 

C 2 (L) 



2e 2 



4Ci(£) 

S 

EVIR2I- 



(A4) 



To determine dynamical properties of the Coulomb 
crystal, we study now the normal vibrational modes, 
which diagonalize the motion Hamiltonian in the har- 
monic approximation. Remember that these ones are 
given by plane-waves of the following form: 



7 5>"" R ' 



e - ipr P r . 



(A5) 



Let us keep in mind that these are not hermitian op- 
erators (Rq)* = R- q . They satisfy the commutation 
relations: [R q ,P£] = iSq iP 5ij, such that the harmonic 
Hamiltonian is finally written like this: 



= E ^^(Rq) 4 (R-q) J + E i P q P 



-q- 



(A6) 

Where the dispersion relation is the Fourier transform 
of the relative potential, see Eq. (5). Allowed phonon 
wavevectors are q = q^bi + q2h 2 , that is, it is a vector 
of the reciprocal lattice (otherwise the plane-waves do 
not satisfy the orthogonality relations on the triangular 



where R r s = R r — R s . Note that we include in i? a h only 
those terms that have the form i/ ah w X Z 2 + X 2 Z 2 . 
The anharmonic coupling constants given by: 

3 e 2 (B° V 



*4 |R° | 4 



3^ 
8 4 



26 hj - 5- 



R! 



r.s) 



R' 



r.s) 



|R?,sl 5 



• (A9) 



Again, we define relative coordinates R" s = R° — R". 

Since we want to induce a quantum gate between ions 
by means of an internal state dependent force in the axial 
(z) direction, the position of the ions will depend, in the 
adiabatic limit of the gate, on the internal state. This 
induces a coupling between qubits and in-plane modes, 
through the dependence of (A8) on R z . To get an explicit 
form for this coupling, let us first introduce the formalism 
for describing a pushing gate between to nearest-neigbors 
in the Coulomb crystal. 



APPENDIX B: PUSHING GATE IN THE AXIAL 
DIRECTION 

Our goal now is to study the performance of the push- 
ing gate under the presence of anharmonic terms. First 
of all, let us collect all the terms that describe our system: 



H(t) = H Q (t)+H ah , 
Ho(t) = H$l + H t (t) 



(Bl) 
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H(j(t) is the harmonic part of the Hamiltonian, and can 
be solved exactly. We include in it the vibrational Hamil- 
tonian: 



A,q 

as well as the force in the axial direction: 



H f = ]T F(t)Z T tf. 



(B2) 



(B3) 



i=i,2 



Under the set of conditions considered in this work, the 
error in the quantum gate is the sum of two independent 
contributions: (i) nonadiabatic effects due to the excita- 
tion of phonons in the axial direction (common to any 
implementation of the pushing gate), and (ii) the errors 
induced by anharmonic couplings. Our strategy will be, 
first, to solve Ho(t), then to study iJ a h in the interaction 
picture with respect to H (t), and finally to evaluate the 
effect of the anharmonic couplings by doing perturbation 
theory. 

For later convenience, we choose the following tempo- 
ral profile for the amplitude of the off-resonant standing 



wave: 



F(t) 2 = F 2 P{t), 



P(t) = e 



-r\t\ 



(B4) 



We define an interaction picture with respect to Hah, 
which has the only unusual feature that H Q (t) is time- 
dependent: 



Uj(t) = ul(t)u(t), 



9tUi(t) 



U^(t)H ah U (t), 



(B5) 



where U(t), Uo(t), are the evolution operators corre- 
sponding to H(t), H (t), respectively, and U\(t) is the 
evolution operator in the interaction picture. 

We consider the limit /3 z Cl, that is, the Coulomb in- 
teraction between ions is small with respect to the trap- 
ping potentials. The problem could be easily extended to 
the general case, but remember that this condition has 
to be fulfilled to avoid crossings between in-plane modes 
and the axial motion of the ions. In this limit, the evo- 
lution operator consists of qubit-qubit couplings, as well 
as the displacement induced by the force (which can be 
considered independently for each ion): 



U (t) = U g (t)e 



-iHi°lt T, j=12 (r,(tra^ :i - V (t)al :/ )c 



u E (t) 



(B6) 



where rj(t) determines the time-dependent displacement 
of each ion's position: 



= H) f 

J — i 



F(t)Z 



r dr, 



(B7) 



Furthermore, we consider the adiabatic limit, r < w z , in 
which the coupling J(t) is given by: 



J{t) = 2/3 z 



ftw z , 



(B8) 



and the ions 1 and 2 are displaced by the state-dependent 
force in the following way: 

Ul{t)ZjU {t) = Z j {t) + 2r ] {t)Z a z j (B9) 
with r)(t) given by: 

F(t)Z 



v(t) = 



(BIO) 



The coupling of the qubits with the in-plane vibrational 
modes is obtained formally by studying the evolution of 
.Hah in the interaction picture: 

Hl h (t) = Ul{t)HMt) « 



E 



- D ri,r 2 l xv ri,r 2 



)*(R-ri,r 



H^(t). (Bll) 



In the last equation we have neglected all the terms that 
depend on the operators Zj , because they give negligible 
contributions in the adiabatic limit of the quantum gate. 
In this limit, Eq. (Bll) reduces to a coupling between 
in-plane vibrational modes and qubits, which does not 
involve any quantum dynamics in the axial motion. Re- 
call that Eq. (Bll) defines H^y C , that is, the coupling 
that we have used in the main text. 

We have already stated the problem in terms that are 
suitable to study the fidelity of the pushing gate. For this 
task, let us assume that the internal states are initially 
in a pure state: 

|*) = Eca|a>, |a) = |00),|01),|10),|ll). (B12) 



In this appendix we are describing both the error induced 
by in-plane modes, and the error due to nonadiabatic 
corrections in the axial motion, so that we study the den- 
sity matrix of the whole system (axial modes, in-plane 
modes, and qubits). Let us consider, an initial density 
matrix for the system, which describes the initial qubit 
pure state, and thermal phonon states p°, p° y : 

Pl = \y)(*\®p° z ®pl y (B13) 

After the action of the quantum gate, the final density 
matrix is: 

ps = E c «4 W(p\ul) {uMuft) (u°p%(u&) 

(B14) 

In Eq. (B14), U", U" y , are the unitary evolution of the 
vibrational modes in directions z, and x— y, respectively, 
proyected into the qubits' quantum state a, that is: 

U« y {t) = e- iH °^T exp^j\a\Hi h (T)\a)dry (B15) 
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where we have used the time-ordered exponential Tcxp, 
to express the evolution operator in the interaction pic- 
ture with respect to i? a h- We define T, the reduced fi- 
delity, as the overlap of the qubits' reduced fidelity ob- 
tained from pf, with a pure qubit-state after the action 
of the unitary operation, U g : 

T=Y J \^\c^T^T^f (B16) 

where we have defined the following quantities: 

= tr.(U?p°(U?)1), 

T^f = tr xy {U? y p° xy (U^) . (B17) 

F z includes the decoherence induced in the quantum gate 
by the finite gate rate, T, something that is general to any 
quantum gate that relies on the adiabatic displacement 
of vibrational modes (see, for example, [12]). Its value 
can be easily evaluated: 

= e -|r,ND| 2 («, + l/2)((K> ct -K>, 1 ) 2 + «<Tj> Q -<<T|>^ 2 )^ 

(B18) 

where t?nd includes the corrections due to nonadiabatic- 
ity, namely: 

fz q r 



Vnd 



-2i 



(B19) 



where Z is the ground state size in the axial trapping 
potential. The contribution to the error can be estimated 
from (B18). For this, we consider that we are in the limit 
of small errors, and minimize T % with respect to all the 
possible initial states, to get the worst-case error, which 
is given by: 



4?7 2 (n + l/2) 



1(( , n a w 



(B20) 



The main point of this work is the calculation of the de- 
coherence that is inherent to this system, that is, the one 
induced by the coupling of the qubits to the in-plane mo- 
tion. Since this one poses a more complicated problem, 
we will be dedicated to it during the next sections to it 



APPENDIX C: DECOHERENCE INDUCED BY 
COUPLING TO IN-PLANE VIBRATIONAL 
MODES 



Note that JF xy is named simply as T in the main text, 
because we focus there on the decoherence induced by 
coupling to x-y vibrational modes only. 

We assume that the ion 1 is at the position R° x = 
(0,0), and 2 is at R°r 2 = (1,0), so that the following 
substitution allows us to express the coupling terms as a 
function of normal modes: 

(R ri y - (R r2 r = \y. «Y ( J - el91 K ( C3 ) 

q,A 



We rewrite here H!^l in a form that is more suitable to 

x y 

study the scaling of the error with the different parame- 
ters of the problem: 



H${t) = FP(r)^F q ^(r) + 



q,k,A,Y 



F$ = e^(l-e^), 



e q , A ei ; A'(l-^)(l-e^).(C4) 



We define: r/(t) 2 = rj 2 P(t), where P(t) is defined in Eq. 
(B4), and r\ = 77(0) = FZ /hu) z . In (C4), we have defined 
the adimcnsional coordinates: 



(C5) 



where we have defined the ground stat e size in the radial 
trapping potential, X = yjh/(2muj r ). The spin-phonon 
coupling constants are: 



G 



;2 

,3£ 
84 



-Zf3 z rfu z 



d 



.(C6) 



Let us see how to express T xy in a form that is suit- 
able for doing perturbation theory in the anharmonic 
couplings. First of all, note that the only matrix ele- 
ments {c^H^a) which are different from zero are those 
with I a) = |01), 1 10) : 



H 



ah 

xy 



(Cl) 



The only quantity that we have to evaluate to calculate 
the contribution to T in (B17) is the following one: 



tr 



xy 



Texp 



H^(r)dr 



(C2) 



The advantage of expression (C6) is that it shows ex- 
plicitly the scaling of the coupling of the qubits with the 
in-plane vibrational modes in terms of the relevant quan- 
tities of the problem. 

In the Taylor expansion of (C2), only time-ordered 
terms appear. Even when we are at finite-temperature, 
Wick's theorem applies [17], something that can be 
shown by means of the path integral formalism, or by 
studying directly the generating functional of the corre- 
lation functions (see the last section of this appendix). 
In this way, one can express every term in perturbation 
theory as integrations of the contractions of the theory, 
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which wc define here in the following way: 

^;- k (i.2) = ^' k (*i-*2). 

K'U^ = (T{^(r)^ k (0)}) 

6{r)e-< T + 9(-r)e^ T ) 5 q ,- k *A,v. (C7) 

where we made explicit that in the limit of high temper- 
atures time-order does not matter and our task should 
be infinitely simplified. 

When calculating terms in perturbation theory, 

1 f°° 

tfqi-qto = ^ J e*»TD%\(u>)du>. (C8) 
We will make use the following identity: 



e^9(±r) = — I . , A , . ; >L: (C9) 



2?r« J-oo w-u^T^ 



So that: 



D x,x 
1,-1 



( W ) 



n A , q 27T(5(w - Wq) + n A , q 27r(5(a) + w*) 



+ 



2^ 



w 2 - (w*) 2 + i£ 



(.CIO) 



The reader is going to find out the reason for our 
choice of exponential pulses: their Fourier transform is a 
Lorentzian, which is better suited for integrations with 
the phonon propagator in Eq. (CIO), 



P{t) = due^Piw), 
2T 



P(u>) 



lo 2 + T 2 ' 



(Cll) 



The energy scales of this problem are J, T, Lu ri and T. 
J and T are related by the requirement that a sign gate 
is realized by the coupling to the axial modes, that is 
J = jT, something that allows us to rewrite the coupling 
constants in a self-consitent way: 



F = 
G = 



3 fX 



do 



(C12) 



Thus, it is clear from the previous expression, that F, G 
are small, both if Xo/do <C 1, and in the case in which 
the quantum gate is also adiabatic with respect to the 
in-plane motion. However, the scaling of the coupling 
constants with these parameters does not allows us to 
extract definitive conclusions about the decoherence, be- 
cause the energies of the vibrational modes are strongly 
corrected by the Coulomb energy (see Fig. 2). 



APPENDIX D: PERTURBATION THEORY 



In this section, we present the calculation of the con- 
tributions to the expansion of JF xy which are of fourth 
order in (X /d Q ), that is, up to C(F 4 ), 0{G 2 ), 0(GF 2 ). 

We have to evaluate time ordered averages of products 
of operators of the form: 



(T{R£ i (t 1 )...Rt n (t n )}). 



(Dl) 



The average () is understood as a thermal average over 
the in-plane vibrational modes. The calculation of time 
ordered correlation functions becomes simplified by the 
use of Wick's theorem. Even when it is usually stated 
in terms of an operator identity that is useful only at 
zero temperature, Wick's theorem can be applied to finite 
temperature, by using the fact that a thermal state is 
gaussian, and thus correlation functions can be obtained 
as functional derivatives of a gaussian functional. We give 
more details of this point in the following subsection. For 
the moment, let us recall that Wick's theorem states that 
(Dl) can be expressed as a sum of products of two-time 
correlation functions, each of them corresponding to a 
possible pairing of the operators: 



(T{Rt 1 (h)...Rt n (t n )}) = 



D 



k„_i,k n (^n-1 t n ) + 
Dtt 3 (tl-t 3 )D^t 4 (t2-U) • • • ifc&ftn-l-tn) + 

. . . (all possible pairings). (D2) 

Note that Wick's theorem allows us to apply the cluster 
theorem to the expansion of the time order exponential 
in Eq. (C2), that is: 



(D3) 



where £^ are disconnected contributions only, that is, 
terms in which the time integrations cannot be factorized 
into separate independent time integrations. 



1. First order in 



The lowest order contribution is given by: 

>nr />oo 

£{1) = H)T2 E G k!- k / dtP{t){Ri{t)R\{t)) 



k,A 



1G 1 

= H)-FT^E G k:"k4k(2nA,k + l). (D4) 



r l 2 



k,A 



where we used x A , k = l/y^w£/wxy, and n A , k is the oc- 
cupation number of the mode with polarization A, and 
wave-vector k. Note that this is purely imaginary, thus 
it should not contribute to the error, after a correction 
of the gate time. 
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2. Second order in H^y 



We find two connected terms: 



q,A 



J P(ti)Dq;\(ti - t 2 )P(t2)dhdt2. (D5) 
We solve the integration in frequency space: 



/ 
/ 



P{t 1 )P{t 2 )D^_ q [t 1 -t2)dt 1 dt2 

(2Tf 



— 2 Iy \. \ 



(2n A , q + 1)^ 



-2 ( U q 1 2tJ q r \ 

_2iXA - q ^Tp + K) 2 + (f 5 + K) 2 ) 2 J ' (D6) 

Since (F/T) fa Xo/do, it is convenient to rewrite this ex- 
pression in the following form. The second disconnected 
contribution of order two in H^y, is: 

q,k,A,A' 

/^P(- + W 0X- q (-)<- A k (-')-(D7) 

The integrations are easily solved in the complex plane by 
contour integration, however, the resulting expressions 
are cumbersome, and bring very little to the discussion 
of the results. The limit of high temperatures, which 
is valid for most of the temperature range considered in 
Fig. 3, is more easily handled: 



£( 3) = --1g 2 £ 



q,k,A,A' 



riA iq nA',k (2P(^ + a/*') + 1P(w\ - a/*')) . (D8) 

Note that f (Wq — uj£ ) contain terms that will diverge in 
the limit r — ► 0, so that, their numerical evaluation is 
necesary for a quantitative understanding of the dissipa- 
tion during the quantum gate. 



3. Third order in 

The only disconnected contribution is given by: 

^ 4) = 4^E f -Wx 



A, A' 
q,k 



Again, we writte here explicitly the result for high tem- 
peratures, that is: 

£(4) = ^F^^^^gW'^^^ x 



A. A' 
q,k 



2P(^)P(^') (p( w * + ug) + P(^ - u£)){mo) 



Where we find, again, a divergence due to resonant terms 
in the limit r — > 0. 



APPENDIX E: SHORT NOTE ABOUT THE 
WICK'S THEOREM AT FINITE TEMPERATURE 

Since it seems that there is some confussion in the lit- 
erature about the validity of this approach, we show here 
that Wick's theorem can be justified at finite tempera- 
tures. 

Let us define the generating functional, which func- 
tional derivation yields the time ordered correlation func- 
tions: 



(T{R kl (t 1 )...R ]i Jt n )}) 
S S 1 



Z\3{x)\, 



Z = tr{e _ ^ Te i £ k ^- 7 ' k(t)Xk(t)dt }. 



(El) 



This could be the starting point for defining a path in- 
tegral representation, but it is much easier to calculate 
explicitly the time-ordered exponential. First of all a 
small detail: it is convenient to consider jk(x) as be- 
longing to the set of continuous functions that satisfy 
j-k(x) = ik(z), so that the exponent is anti-hermitian. 
Note that the generating functional (El) has the same 
form as the evolution operator of a forced harmonic os- 
cillator, thus, its form can be explicitly obtained, by fol- 
lowing the same lines. In any case, since the conmutators 
of the exponents at different times, are itself quadratic 
forms of the currents j'k, it is clear that the generating 
functional is gaussian. Since its second functional deriva- 
tive yields the time ordered correlation function of tow 
position operators, the only choice is the following one: 

Z\j{x)]/Z, = (E2) 
cxpj^ (+i J j_k(ti)Dk(*i - *2)jk(*2) dtxdt^j 



where the kernel of the functional, Dk{t\ — t 2 ), is the 
correlation function defined in the previous sections. 
Higher order functional differentiation leads immediately 
to Wick's theorem. 
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APPENDIX F: EFFECTIVE SPIN-SPIN 
INTERACTIONS INDUCED BY A 
WALKING-WAVE 

Finally, we introduce here a scheme for inducing ef- 
fective spin-spin interactions which relies on the use of 
a walking wave. This is a different approach to the one 
presented in [19], which relies on the coupling of the in- 
ternal states of the ions to the motion by means of a 
standing wave. The latter approach might be difficult 
to implement, since one needs to fix the position of the 
ions relative to the minima of the standing wave. This 
condition is not necesary in the walking wave scheme, 
something that may be useful in quantum simulations 
with ions in Paul or Penning traps. 

A walking wave, like the one used by the NIST group 
for the geometrical phase gate [3] , leads to a spin-motion 
coupling of the following form: 



n ( 

2 l f 



ik(xQ-\-Xj) — iujL t 



e -ik(x -\-Xj)+iuj L t . p. 



(Fl) 



where u>l is the detuning of the laser beams that create 
the walking wave. We consider that x is a direction trans- 
verse to the Coulomb crystal, that is, the radial direction 
in the case of a linear chain, or the axial direction in a 
2D Coulomb crystal. If the walking wave propagates in 
this direction, then all the ions are at the same position, 
xq with respect to the walking wave. Note also that we 
are neglecting ac Stark shifts, which can be cancelled in 
a laser configuration like the one used at NIST [3]. 

We consider here that the vibrational modes of the ions 
are in the "stiff" limit, that is, the trapping potential is 
much larger the Coulomb repulsion between ions. The 
validity of this limit can be quantified by means of the 
following quantity: 



0x 



(F2) 



where do is the mean distance between ions. If j3 x <C 1, 
then the energy of the vibrational modes is not strongly 
modified by the Coulomb interaction, and all the modes 
are close in energy. Indeed, (3 X is an approximation to 
the dispersion in the energies of the vibrational modes, 
lu u . Let us say that uj x is the trapping frequency. Then, 
all the energies lie below u) x , and satisfy that: 



(1 - Px)U) x < LU n < LU X 



(F3) 



If we consider the Lamb-Dicke expansion of (Fl), then: 



ikxQ — iujLt 



+ h.c.) a z 



ike lkxo e~ lULt J2 x j + h - c - a ]- ( F4 ) 



Our idea is the following: in the case (3 X <C 1, one can 
tune to the red-sideband with respect to all the vi- 
brational frequencies in (F4) (see Fig. 4), because in this 



D>0.6 
<D 
c 
HI 



CO, 



10 15 

mode number 



FIG. 4: Radial vibrational normal modes of a chain with 
20 ions, and ration between radial and trapping frequencies: 

U) r /u)z = 10. 



case, ions are close to be independent. In the end, this 
will lead to a coupling of the form J2j x j a ji which will 
induce an effective spin-spin interaction just in the same 
way as in our previous proposals, with the advantage that 
the amplitude of the force does not depend now on the 
position of the ions relative to the phase of the walking 
wave. 

Let us describe this situation formally. First, remem- 
ber that the position of the ions can be expressed in terms 
of normal modes: 



1/2 



(F5) 



where M. are matrices that diagonalize the harmonic mo- 
tion of the ions (they are basically sines and cosines) . 

We want to tune u>l to the red sideband with respect 
to all the modes u> n , so that we choose lul < w n , and 
uj x — u> n -C ui n . Under these conditions we can keep only 
the red sideband terms of (F4): 



9. 



ikxQ 



2mu! 7 



1/2 



Mj, n a*. (F6) 



The rotation of the terms in (F6) can be eliminated 
by considering the phonons in a rotating frame, in which 
their energy is shifted by u>l. The whole hamiltonian in 
this frame reads: 

H = H o + H f (0) 

= ^ hS nala n + ^ {9n a n + 9nal) , (F7) 
n n 

with S n — oj n — ujl- Note that (F7) induces a state de- 
pendent force on all the ions simultaneously. A difference 
with the pushing gate approach is that vibrational ener- 
gies suffer a shift that we can controll by tuning uil- 
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To deal with this coupling we use, as usual, a canonical 
transformation, to map our system to a quantum spin 
model: 



E 



{if n a n - rj„al) 



9n 



(F8) 



Under this transformation, the internal states of the ions 
follow the dynamics of a spin-spin Hamiltonian: 



H - e^He* = H + -J . k ^ 0*0%. 



The effective spin-spin interaction is given by: 



Jj.k = > Z~Mj n Mk,n- 



(F9) 



(F10) 



where is the force induced by the walking wave, F — 
ik(Sl/2)e~ lkx " . The spin-spin interaction does not de- 
pend on the absolute position of the ion crystal, or the 
ion plane. 

Note that one important difference with our previous 
work is that now the shifted frequencies appear in the 
denominator, thus, we can expect that the dependence 
of Jj.k on the inter-ion distance will be different, and 
even suitable to be controlled by tuning W£. 

Let us now estimate the strength of Jj t k in the simplest 
case. For this, let us recall that the normal vibrational 
modes can be written like this: 



K) 2 = (w x f(l + p x V n ), 



(Fll) 



where V n is the contribution from the Coulomb interac- 
tion. So to say, V n contains the dispersion of the modes 
in Fig. 4. V n can be obtained from the harmonic terms of 
the Coulomb interaction by means of the transformation 
defined by the matrices M: 



M jtn V n Mk,n 



1 



(F12) 



On the other hand, the shifted modes fulfill that: 

5 n = (u x )(l + p x V n ) 1/2 -u L . (F13) 
If these two conditions arc met: 



Px < 1, 
P x u x /(uj x - uji) < 1, 



(F14) 
(F15) 



then we can express u> n and S n in a series in f3 x , and 
get the following approximate relation for the effective 
interaction: 



Jj,k ~ Px 



= Px 



1 



2m(uj x - u L ) 2 \j - k\ 3 



Fx 

x -u) L J \j-k\ : 



UJ X 



(F16) 



where we have also used that lu x — lul <C u> x , and xq is 
the ground state size in the trapping potential. Note that 
condition p x <C 1 has to hold if we want to be able to 
resolve the red sideband simultaneously for all the vibra- 
tional modes, but condition (F15), on the contrary, can 
be avoided, the only difference being that the effective 
spin-spin interaction will deviate from the 1/r 3 decay in 
Eq. (F16). 

The expression (F10) does not tell us much about the 
interaction rate that can be achieved in an experiment. 
For this, one has to take into account the error that is 
introduced in the simulation by the coupling of the qubits 
(effective spins) to the vibrational modes. 

We have already studied in detail this error (see [19]), 
which depends on the displacement of the modes. The 
error in the quantum simulation is induced by the entan- 
glement of the vibrational modes with the effective spins, 
and it can be approximated by the following quantity: 



V = 



\Fx f 



h(uj x - uj l ) ' 



(F17) 



which is simply an approximation to the coefficients rj n 
in Eq. (F8). The strength of the effective interaction is 
J w f3 x ri 2 uj x . Thus, we find the same relation that in the 
case of the standing wave. 

Last, one has to take into account the effect of the 
first term in Eq. (F4). Even in the limit fi/w£ <C 1, 
one could get into troubles in case that corrections of the 
form f2 2 /u>l appear. We will see that this is not the case, 
because this term is diagonal in a z . It can be rewritten 
in the following way: 



H c = e(t)a z , 
e(t) = Q cos(kx — u>Lt), 



(F18) 



such that the eigenstates are the same, with time evolu- 
tion given by: 



!*(*)> 
c T (t) 

c T (t) 



c T (*)l T>+ci(t)| I), 

e -i J* e(t')dt' ^ 

g— i(n/wi,)(sin(/cxQ— wi,t)— sin(fexo)) (F19) 



Thus, the effect of H c affects terms of the form Ba x , 
and operators a\ get corrections of order BSI/ujl- This 
means, that condition O/w; <C 1 is enough to guarantee 
that these corrections are small compared to the effective 
spin-spin Hamiltonian interaction strengths. 
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